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ABSTRACT 

We use detailed Lya radiative transfer calculations to further test the claim of Rauch 
et al. (2008) that they have detected spatially extended faint Lya emission from the 
elusive host population of Damped Lya Absorption systems (DLAs) in their recent 
ultra-deep spectroscopic survey. We investigate the spatial and spectral distribution 
of Lya emission due to star-formation at the centre of DLAs, and its dependence on 
the spatial and velocity structure of the gas. Our model simultaneously reproduces the 
observed properties of DLAs and the faint Lya emitters, including the velocity width 
and column density distribution of DLAs and the large spatial extent of the emission 
of the faint emitters. Our modelling confirms previous suggestions that DLAs are 
predominately hosted by Dark Matter (DM) halos in the mass range 10 9 5 — 10 12 M Q , 
and are thus of significantly lower mass than those inferred for L„ Lyman Break 
Galaxies (LBGs). Our modelling suggests that DM halos hosting DLAs retain up to 
20% of the cosmic baryon fraction in the form of neutral hydrogen, and that star 
formation at the centre of the halos is responsible for the faint Lya emission. The 
scattering of a significant fraction of the Lya emission to the observed radii, which 
can be as large as 50kpc or more, requires the amplitude of the bulk motions of the 
gas at the centre of the halos to be moderate. The observed space density and size 
distribution of the emitters together with the incidence rate of DLAs suggests that 
the Lya emission due to star formation has a duty cycle of ~ 25%. 
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1 INTRODUCTION model for structure formation. The model assumes a simple 

i— : ri i -i _ I i relation between the size of the damped absorption and Lya 

Rauch ct al. 2008, hereafter R08) recently reported the re- . . . , , T , . ., . 

1 J emission regions, and proposes that the L ya luminosity is 

suits of an ultra-deep spectroscopic survey for low surface ,. , , , , , , , , , InTTnrJ c ^ u 

, , r « , i proportional to the total halo mass. BH09 further corrobo- 

brightncss Lya emitters at redshilt z ~ 3. A 92 hour fong ex- , , x , ,. ,. . , , 

. , . , , , rated the suggestion that cooling radiation is not expected 

posure with the LbO VLi tORbi instrument yielded a sam- , , . . .„ , , JT . . , 

J _ 1S to contribute sigmhcantly to the observed Lya emission, and 

pie of 27 faint line emitters with fluxes of a few times 10 , , , ... j_i.ii i i_ e 

_2 i i i that the emitters are most fikely powered by star formation, 

erg s cm , which they argue are likely to be dominated by T , , , , „ T . , . , . , , „, , . , 

fn the modef, DLAs are small galaxies hosted by DM halos 

Lya. They further conclude that the large comovmg number .., . lr ,9 5i 1^12 »,r ji ., 

J J _,, _., with masses m the range 10 to 10 Mq and have rather 

density, 3 x 10 n.™ Mpc , and the large covering factor , , r , . , , T ,, ,., , , ,.. 

. r ' a large low surface brightness Lya hafos which extend to radii 
dN/dz ~ 0.2 — 1 suggest that the emitters can be iden- 
tified with the elusive host population of damped Lya sys 



of up to 50kpc or larger. 



terns (DLAs) and high col umn densi ty Lym an limit systems. In order t o fit the observed size distribution of the faint 

iBarnes fc HaehneltJ (|2009l . h ereafter |BH09T|. building on the Lya emitters, |BH09| assumed that the Lya emission ex- 
successful m odel for DLAs of lHaehnelt. Steinmetz. fc Rauchl tends to radii somewhat larger than is required to repro- 
(| 19981 . bOOd ). presented a simple model that simultaneously duce th e incidence rate for DLAs. However, no modelling of 
accounts for the kinematic properties and incidence rate of the gas distribution and Lya radiative transfer was done, 
the observed DLAs and the luminosity function and the size We present such modelling here to investigate whether the 
distribution of theHol emitters in the context of the ACDM sizes, surface brightness profiles and spectral line shapes can 

be reproduced with simple but plausible assumptions for the 
distribution and the physical properties of the gas in the DM 
* E-mail: lab@ast.cam.ac.uk (LAB); haehnelt@ast.cam.ac.uk halos suggested by our previous modelling to be the hosts 

(MGH) of the faint emitters. 
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Studies of the radiative transfer of Lya have a long 
history. The first investigations employed app roximate 
calcu l ations and simple physical arguments jZanstral 
1 19491 ; lUnnol Il952l ; lOsterbrockl Il962l ; lAdamsl Il97ll . Il972l ). 
Later, analytic solutions were found for sim ple geome- 
tries i n the l i mit o f large optical depth (jHarringtonl 
1 19731 ; iNeufeldl Il990h . It was realised early, however, 
that Monte Carlo simulations provided the most flexible 
method for investigating arbitrary geometry, density 
distributions and velocity structures. Many investiga- 



1972) 



1979) 



tions hav e em plo yed s uch techniques - | Averv fc H ouse 
' 19681 ); lAueri lfl96sl); ICaroff. Noerdlinger, fc Scargle 



Panagia fc Ranieri 



Natta & Bcckwith 




IE. 
1 



Bonilha et al 



Ahn. Lee, fc Lee 



(2000 . l200ll. 120021); IZheng & Mira lda-Escudel (12 002): 



(2006, 



Tasitsiornil 



here- 



ICantalupo et all (120 05); Han sen fc Ob] 
(2006), Diikstra, Haiman, & Spaans 
after |PHS06h. IVerhamme. Schaerer fc Maselhl l|2006l ); 
lLaursen. Razoumov. fc Sommer-Larsenl l|2009l ). 

Our modelling is most sim i lar to that of DHS06 and 
IVerhamme. Schaerer. fc Masellil {2006), who modelled Lya 
radiative transfer in collapsing proto-galaxies and high- 
redshift galaxies respectively. In many instances we make 
similar assumptions as these authors, but our code was de- 
veloped independently. 

The paper is structured as follows. In Section 2 we give 
a brief summary of the salient features of our Monte-Carlo 
code for Lya radiative transfer and show results for standard 
test problems. In Section 3 we discuss our assumptions for 
the distribution and the physical properties of the gas. We 
also show the dependence of the surface brightness profile 
and the spectral line shapes on these assumptions. In Section 
4 we present the results for a consistent model of the size 
distribution and the luminosity function of the faint Lya 
emitters. The technical details of the implementation of Lya 
radiative transfer Monte-Carlo algorithm are described in an 
appendix. 



2 THE MONTE-CARLO RADIATIVE 
TRANSFER CODE 

2.1 General properties of the code 

Lya is a resonant transition in hydrogen. The optical depth 
at line centre to of a typical Hi region is thus generally very 
large. Lya photons, which are expected to be produced copi- 
ously in star- forming galaxies, will therefore undergo many 
scatterings. Escape from the Hi region usually requires dif- 
fusion in both frequency and space. 

The properties of the emergent radiation depend sensi- 
tively on not only the spatial distribution of the gas, but also 
on its kinematics, ionization state, temperature and dust 
content. With modern computers Monte-Carlo sampling of 
the diffusion process has become the method of choice. We 
employ such a code for simple spatial and kinematical con- 
figurations. 

Our code makes a number of well-tested assumptions 
and approximations. Throughout, we will use a dimension- 
less frequency variable x, the frequency displacement from 
line centre in units of the Doppler frequency width (see ap- 
pendix A for details). As the photons experience a large 



number of scatterings before they escape, we can inject all 
our photons at line centre in the fluid frame of the gas 
(^initial = 0) for simplicity. The photon scattering is par- 
tially coherent — in the rest frame of the scattering atom, 
the final frequency of the photon differs from the initial fre- 
quency only by the (often negligible) effect of atomic recoil. 
The direction of the photon post-scattering is chosen from 
a dipole distribution, although choosing an isotropic distri- 
bution makes little difference. We also incorporate a cosmic 
abundance of deuterium (see |PHS06l) . Finally, we give all 
results in the rest frame of the centre of mass of the emitter. 

Monte Carlo Lya RT codes like the one used here can be 
significantly accel erated by skipping scatt er ings in the core 
of the line-profile |Ahn. Lee fc Leei (|2002l ). |PHS06| . In the 
Doppler core, the optical depth is so large that the photon 
experiences very little spatial diffusion while it awaits an 
encounter with a rapidly moving atom that will scatter it in 
frequency space into the wings of the line profile. The gas 
distributions we consider here are spherically symmetric. We 
exploit this symmetry by using an arrangement of spherical 
shells. 

A detailed description of the implementation and rele- 
vant formulae of our code can be found in Appendix lAl 



2.2 Testing the code 

There are a number of analytical solutions to Lya radiative 
tra nsfer probl e ms ag ainst which we have tested our code. 

iHummen (|1962l ) calculated the redistribution function 
for the case of coherent scattering in the atom's frame with 
radiation damping (i.e. incorporating the L orentzian natur al 
line width). The result (Equation 3.12.2 of lHummeri (|l962l )) 
is Ru-b{x, x') Ax' , which is defined as the probability that 
a photon, whose frequency prior to absorption was x, is re- 
emitted with frequency in the range [x' ,x' + dx']. In the top 
left panel of Figure [T] we compare the analytic formula with 
the output of our code. The thin, smooth line is the analytic 
solution and the thick histogram is the output from our code 
for x — 0, 2, 5 (blue, red, black respectively) for T = 10 K 
(a = 0.0149). The two lines are almost ind istinguishable. 

iHarringtonl (|l973l ) and INeufeldl (|l990h derived an ana- 
lytic expression for the spectrum J(x, to) of radiation emerg- 
ing from an optically thick (^/ttto > 10 3 /a), uniform, static 
slab of neutral hydrogen, where line-centre photons are in- 
jected at the centre of the slab, atomic recoil is neglected 
and to is the line-centre, centre-to-edge optical depth. The 
solution normalised to (47r) _1 is given by, 



J(a;,T ) = 



V6 



240FaT O coshjy 

7r 3 /54 |a; 3 |/aTo 



(1) 



The comparison of our code with this analytic solution is 
shown for T = 10K in the top right panel of Figure [1] The 
agreement is again very good, especially as to increases. 

Next we compare the results of our code to the analogue 
of the above solution for a uniform, static sphere (|PHS06l . 
Equation (9)). As the bottom left panel of Figure Q] shows, 
the agreement is also very good. 

To test the code for the case of bulk m otions in the gas, 
we com pare our code with the modelling of lLoeb fc Rvbickil 
(1999). Loeb & Rybicki investigated Lya RT for a uni- 
form Hi distribution undergoing Hubble expansion with 
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Figure 1. Four tests of our Lya RT code. Top left panel: The redistribution function R\i_^(x,x') gives the probability that a photon, 
whose frequency prior to absorption was x, is re-emitted with frequency in the range [x' , x' +dx']. The thin, smooth curve is the analytic 
solution and the thick histogram is the output from our code for x = 0, 2, 5 (blue, red, black respectively), and T = 10 K (a = 0.0149). 
The analytic solution and our Monte-Carlo results are almost indistinguishable. Top right panel: For an optically thick (y/nro > 10 3 /a), 
uniform, static sla b of neutral hydro gen, where line-cent re photons are injected at the centre of the slab, we can compare with the 
analytic solution of iHarringtonl lll973t) and iNeufeldl dl99(Jl . We set T = 10 K, and to is the line-centre, centre-to-edge optical depth as 
labeled on the plot. The solid histogram from our code agrees again well with the analytical solution, especially as to increases. Bottom 
left panel: As above, but for a uniform sphere. The agreement is very good. Bottom right p anel: Ly a RT through uniform Hi undergoing 
Hubble exp ansion with velo c ity vp = H(z)r. The histogram is the result from the code of IDH S06. the green solid curve is the solution 
of by LR99 lLoeb "fc Rvbickil i 19991 ), and the dashed blue line is the result of our code. The agreement is again excellent. 



Vb = H(z)r, where H(z) is the Hubble constant at redshift 
2. For photons blueward of line-centre in the fluid frame, the 
expansion of the universe will eventually redshift the pho- 
ton back into resonance. Only photons redward of line centre 
can propagate to infinity. The optical depth to infinity is 
Tim (To a Vth _ x„ 



(2) 



s/ttH{z)x x 
where x is the frequency in the fluid frame. This defines a 
critical frequency x t . Photons with x -C x„ are redshifted 
enoug h to stream freely. We follow here lLoeb fe Rvbickil 
l|l999l) and use a new frequenc y variabl e x = x/x*. We refer 
the reader to Appendix B2 of DHS06 for a careful discus- 
sion of the modifications necessary to permit a meaningful 
comparison between the two codes. 



In the bottom right panel of Figure [T ]we compare th e 
results of our code with those of lLoeb fc Rvbickil l|l999T i. 
and[DHS06. The relevant parameters are: z = 10, nm = 2x 
10 _7 (1 + z) 3 cm" 3 and T = 10K. Note that / P(x) dx = 1. 
The agreement is again excellent. 
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3 Lya RADIATIVE TRANSFER IN DLAS / 
FAINT Lya EMITTERS 

3.1 A simple spherically symmetric model for the 
spatial distribution and kinematics of the gas 

3.1.1 Observational and theoretical constraints 

Our knowledge of the spatial distribution and kinematics 
of the neutral hydrogen in the DLAs / faint Lya emitters 
is still somewhat limited. The statistics of the occurrence 
DLAs and their column density distribution give us inte- 
gral constraints on the spatial distribution of the gas, while 
the velocity distribution of low ionization species tracing 
neutral hydrogen gives us some indication of the velocity 
range of bulk motions. The bulk motions of the gas ap- 
pear to have velocities tha t range from a few tens of km/s 
to several hundred km/s l|Prochaska fc Wolfel [1997V The 
relative contribution of ordered and random motions and 
the role of gas inflow and galactic winds is, however, still 
very uncertain (Haehnelt, Stcinmctz, & Ranch 1998, 2000; 



Razo umov et al 
20091 ). 



20081 ; iPontzen etail 120081 ; ITescari et ail 



If the identification of the faint Lya emitters as DLA 
host galaxies is indeed correct, then this gives us for the 
first time constraints on the spatial extent of the gas dis- 
tribution in individual objects for a sizable (Lya) emission 
selected sample. While there is still significant ambiguity 
due to th e unkn own "duty cycle" for Lya emission (as dis- 
cussed by BH09) , theoretical modelling gives us a handle on 
the halo masses and virial velocities expected to host the 
DLAs / faint Lya emitters. By duty cycle, we mean the 
fraction of DM host haloes detectable through their Lya 
emission at any given time. There are a number of reasons 
why this fraction may be smaller than unity. The haloes, for 
example, may not continuously form stars, but geometrical 
factors due to a preferential escape of the Lya emission in 
certain directions due to variations in the neutral hydrogen 
column dens ity (and/or the dust column density) shou l d also 
contribute llLaursen. Razoumov. fc Sommer-Larsenl 120091 : 
ILaursen, Sommer-Larsen. fc Andersenll2009f r 

Giv en these uncertainties, we have decided to follow 
DHS06 and investigate Lya radiative transfer for sources 
at the centre of DM halos with a range of masses, a duty 
cycle for Lya emission, a simple spherically symmetric gas 
distributions without dust, and either inflow or outflow with 
velocities which vary as a power law with radius. This al- 
ready leads to a rich variety in the predicted surface bright- 
ness profiles and spectral shapes and allows us to study the 
influence of important physical parameters. 



3.1.2 The assumed radial distribution of neutral hydrogen 

We begin by specifying the radial distribution of neutral 
hydrogen in a given halo. The total amount of hydrogen is 
set relative to the cosmicQ mass fraction of hydrogen /h = 
On/flm- Throughout we assume a helium fraction of Y p — 
0.24. There are a number of reasons to suspect that the 
hydrogen mass fraction in a typical halo is lower than /h- 



1 The relevant cosmological parameters used in this work are: 
(h,n M ,n b ,n A ,as,n) = (0.7,0.3,0.045,0.7,0.9, 1). 



Firstly, baryons are subject to the smoothing effects of gas 
pressure. Secondly, gas that forms stars is both ionised and 
extremely compact. Stars will also ionise the neutral gas 
around them — this is a source of Lya but also reduces 
the amount of Hi that remains to scatter photons. Finally, 
stellar and AGN driven galactic winds are expected to drive 
gas out of galaxies back into the IGM. As a first attempt at 
modelling this effect, we reduce the total mass of baryons in 
the halo to a fraction f c of the cosmic value. 

On top of reducing the amount of neutral hydrogen in 
a typical halo, it is known that the UV background will sig- 
nificantly ionise the gas in halos too small to self-shield. It 
is also easier for galact ic winds to drive gas out of small, 
shallow halos. In lBH09l . this effect was implemented via an 
exponential suppression of the cross-section of neutral hy- 
drogen below a critical circular velocity u c ,o. This was nec- 
essary in order to fit the observed velocity width distribution 
of associated low ionization metal absorption. Here, we will 
implement this suppression by reducing the total amount 
of neutral hydrogen in halos below v c ,o, such that the total 
mass of neutral hydrogen in a halo of mass M v is 



Mhi = /e/H exp 



Vc,0 
Vc 



(3) 



(see IBH09I for a detailed discussion). We use the fiducial 
parameters a c = 3, v c> o = 50 km s . 

For the radial distribution of the gas, we assume an 
NFW pro file at z = 3, characterise d by a scale radius r a 
defined in iNavarro. Frenk. fc Whitel ll 19961) . Following the 
simulations of Mailer fc Bullock! 1 2004 Equation (9)), we 
alter the NFW profile to give the halo a thermal core at 
~ 3r s /4. The profile is then specified by the total mass of 
the halo M v and the concentration parameter c v = r v /r s . 
For dependence of the concentration parameter on the mass, 
we take the mean valu e of the c v — M v correlation as given 
by IMaccio et al.l (|2007l ). 



C v = Co 



10 n M Q 



1 + z 



(4) 



For dark matter, IMaccio et al.l |2007l ) found that cq w 3.5, 
with a log-normal distribution and a scatter around this 
m ean value of A ( lnc v ) = 0. 33, in agreemen t of th e results 
of iBullock et al.l l|200ll ) and IWechsler et alj (|2002l ). As we 
will find later, a significantly larger Co is appropriate for 
the baryons; we will use the column density distribution of 
DLAs to constrain Co in Section [4.21 As we discuss further 
in Section 14.51 the gas in the DLAs / faint emitters can 
be expected to self-shield against the meta-galactic ionizing 
UV background at z ~ 3. The corresponding self-shielding 
radius in the DM halos we are studying here is generally 
smaller than the virial radius. We therefore set the outer 
radius of the Hi to be the virial radius in our modelling 
and also ignore radiative transfer through the IGM. We set 
T = 10 4 K as a fiducial temperature. 



3.1.3 The assumed kinematics of neutral hydrogen 

The biggest uncertainty in our modell ing is p robably the 
kinematical state of the gas. We follow DHS06 and for our 
fiducial model we assume the gas to be infalling with a 
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power- law radial velocity profil^l parameterised by u amp and 



impact parameter (") 
2 3 4 



Ubulk(r) = -Damp 



(5) 



where r v i s the virial radius , for w hich we follow the def- 
inition of iMaller fc Bullockl { 2004). Values in the range 
a 6 [—0.5, 1] should be reasonable. The upper limit de- 
scribes a spherical top-hat collapse, while the lower limit 
represents the accretion of massless shells onto a point mass 
(v 2 ~ GM/r). A spherical top-hat would have u amp = v c . 

While gas infall will certainly be an important feature 
of the gas kinematics, st ar formation driven outflows are also 
likely to play a role (see lVeilleux. Cecil, fc Bland-Hawthornl 
2005, for a review). Note that for a spherically symmetric 
gas distribution, the red and blue side of the line profiles will 
just be interchanged if the gas is assumed to be outflowing 
instead of inflowing with the same velocity profile (ignor- 
ing the very slight effects of recoil and deuterium). For the 
more massive and more actively star-forming LBGs, galactic 
winds have been su g gested to sweep up an expanding shell 
ilPettini et alj |2000|. |2002|; IVerhamme! Schaerer. fc Masellil 



20061; ISchaerer fc Verhammel 120081 : IVerhamme et all I200I 



Quider et al.l 



2009). 



We investigate such a configuration in 



Section 13.31 



3.1.4 The fiducial model 

IROI and [BH09 identified star-f ormation as th e most likely 
source for the faint emitters. IWolfe fc Che"n1 (I2006T ) used 
continuum emission to place stringent limits on extended 
star formation in DLAs. We have thus assumed a centrally- 
peaked emissivity — all photons are created at r = 0. 

Before we attempt to model the data of iROSl in detail, 
we will consider the effect on the spectra and surface bright- 
ness distribution of altering the parameters of our model. 
Our fiducial model parameters are (z, M v , cq, / c , v amp , T) = 
(3, 10 n M Q , 25.3, 0.2, v c , 10 4 K). The values of c and f c that 
we have chosen will be justified in Section 14.21 The surface 
brightness S scales with the total luminosity Z/Lya, which 
will be given in the caption to each figure. The values of 
the luminosity chosen will be justified in section 14.31 The 
velocity and column density profiles for the fiducial model 
are shown in Figure [21 where the column density is as seen 
along a sightline that passes all the way through the halo at 
an impact parameter y. 



3.2 Lya radiative transfer in individual halos with 
gas infall 

In this section, we will consider the effects of changing the 
most important parameters in the model: mass, concentra- 
tion, velocity profile (parameterized by the power law index 
a) and baryon fraction. Further parameters (temperature 
and amplitude of the bulk motions of the gas) are discussed 
in Appendix [Bl 



2 W e use the modification to this law for a < given in Equation 
10 of DHSol 
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Figure 2. Top: The velocity profile for the fiducial halo, for dif- 
fering values of a as shown in the legend. Bottom: The column 
density profiles for the fiducial model with the total mass as given 
in the legend in units of Mq. The column density is calculated 
along a sightline that passes all the way through the halo at a 
distance r from its centre. Note that the 10 9 Mq line has been 
boosted by a factor of 10 4 to make it visible on the given axes. 
The horizontal dashed line indicates the minimum column density 
of a DLA. 



3.2.1 The effect of halo mass and concentration 

The models of lBH09l give us a handle on the masses of the 
halos that host DLAs. Figure [3] considers halos with masses 
of M v = I0 9 , I0 10 , 10 11 , 10 12 M Q , which at z = 3 correspond 
to v c = 23, 50, 107, 231 km s _1 , r v = 8, 17, 37, 81 kpc, 6 max = 
1.05, 2.25, 4.9, 10.5 arcsec, where fo max is the angular radius 
corresponding to r\. 

We see that, as the mass increases, the emerging pho- 
tons emerge bluer. This is because, as we add more gas, the 
central column densities increase and the photons must shift 
further from line centre in order to escape. The a = 1 pro- 
files are increasingly double-peaked for lower masses, while 
the a = —0.5 profiles only ever have one, blue peak. This 
is because the innermost region of the a = 1 halo has the 
smallest bulk velocity, and thus most resembles the uniform, 
static sphere of Figure [1] As the bulk velocities increase at 
the centre of the halo, the amount of energy transferred be- 
tween the gas and photons in each scattering is increased, 
favouring one of the two peaks (the blue/red peak for in- 
flowing/outflowing gas respectively). 

The surface brightness profile shows that the dominant 
effect in increasing the mass is that the virial radius (which 
we have assumed to be the outer radius of the Hi) increases. 
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Figure 3. Spectra (left panels) 
and surface brightness profiles 
(right panels) for emission in halos 
with total mass as given in the 
legend in units of Mq. The top 
panels are for a = —0.5, with a 
total Lya luminosity which scales 
with the mass of the DM halo 
as required to fit the observed 
luminosity function in Fig. 10 
(see section 4 for more details). In 
the order given in the legend the 
luminosities are L^ ya = (9.6 X 
10- 6 ,0.88,22,240)xl0 42 erg s" 1 . 
The bottom panels are 
for a = 1, with Z/Lya 

. x lO -7 , 0.044, 1.1, 12.1) x 
10 42 ergs -1 . Note that the 
surface brightness for the 
M = 1O 9 M model has been 
raised by a factor 10 4 to be able 
to plot it in the same plotting 
window. As the mass increases, 
the emerging photons emerge 
bluer, and are scattered to larger 
radii in the larger halos. The 
dashed horizontal line is the 
detection threshold of the Rauch 
et al. emitters. 



The a — —0.5 profile is much more centrally peaked than the 
a — 1 profile. The reason is that the larger bulk velocities 
at the centre of the a = —0.5 halo can shift the photon into 
the wings of the spectral line, resulting in reduced spatial 
diffusion. 

For the concentration parameter, we have considered 
the values 1.8,3.5,9.4 and 25.3. The spectra and surface 
brightness profiles for these models are shown in Figure |U 

As the concentration increases, the photons generally 
emerge bluer. The a — 1 profile becomes more double- 
peaked as the concentration increases, because the Hi col- 
umn density increases at smaller radii, where the bulk ve- 
locity is lower. The surface brightness profile is also more 
centrally peaked for higher concentrations, as more scatter- 
ings occur at smaller radii. 



3.2.2 The effect of baryonic fraction/ column density 

The effect of changing the baryonic fraction f c is shown in 
Figure [5] As baryons are removed from the halo, the spec- 
trum shifts toward x — and the surface brightness profile 
becomes more centrally peaked, as there is less gas in the 
outer parts of the halo to scatter the photons. The baryon 
fraction where the photons are not scattered efficiently any- 
more to the virial radius (where we have the gas distribution 
assumed to cut-off) corresponds to a HI column densitj0 of 
about 10 16 cm" 2 . For very small / e and a = —0.5, some of 
the photons can escape the halo without scattering at all, 



3 Which is somewhat dependent on the spatial profile of neutral 
hydrogen and the velocity field. 



creating a very narrow peak at x = 0. Some of the spec- 
tra also show a trough at x ~ 6 due to deuterium. Note, 
however, that for f c <C 0.1 the gas in the halos would not 
be able to self-shield anymore against the meta-galactic UV 
background at z ~ 3 so the model would be internally in- 
consistent for such small baryonic fractions. 

3.3 Expanding Shells 

For the much brighter Lya emission seen in many LBGs, the 
Lya emission is systematically red-shifted by several hun- 
dred km s _1 . This is generally attributed to backscattering 
of the Lya radiation from a w i nd-dr iven expanding shell. 
IVerhamme. Schaerer. fc Masellil (|2006h used a Lya RT code 
similar to ours to calculate the emergent spectrum from an 
expanding shell. We will investigate here in a similar fash- 
ion if this picture could apply to the R08 emitters. We con- 
sider a shell of Hi, where the inner radius is a fraction f r 
times the outer radius r max . The shell is expanding at a uni- 
form velocity V exp , and has a column density of Nm- The 
temperature is set by the Doppler velocitjQ br- Our fidu- 
cial model has the parameters: (r ma x, fr, V e xp, Nm, br) = 
(30kpc, 0.9, 300 km s" 1 ^ x 10 20 cm -2 , 40 km s" 1 ), simi- 
larly to the fiducial model of Verhamme et al. For the outer 
radius of the shell r max we have a chosen a value similar to 
the models of the previous sections. 

Figure|6]shows the effect of altering the parameters that 
Verhamme et al. found to be most important, the column 

4 This is defined in the same way as the thermal velocity dis- 
persion (in the Appendix), with a possible contribution from a 
turbulent velocity dispersion, added in quadrature. 
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Figure 4. Spectra and surface 
brightness profiles for emission in 
halos with concentration param- 
eter co as given in the legend. 
The top panels are for a = 
—0.5, for which L^ ya = 2.2 X 
10 43 erg s — 1 . The bottom pan- 
els are for a = 1, for which 
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Figure 6. Spectra and surface 
brightness profiles for the modelling 
of expanding shells. The panels show 
the effect of altering the parameters 
that Verhamme et al. found to be 
most important, the column density 
and the expansion velocity, as given 
in the legend. Note that, for small 
column densities (in the upper pan- 
els), a certain fraction of photons 
were able to escape directly with- 
out scattering. These photons form 
a delta function at x = 0, and at 
small impact parameter, which has 
been removed from these plots for 
clarity. As the column density in- 
creases (top panels), the spectrum 
becomes redder, and the small blue 
peak disappears. The most signifi- 
cant change in the surface brightness 
profile comes at small column den- 
sities as photons are able to escape 
directly, leaving less to be scattered 
at large radii. As the expansion ve- 
locity increases (bottom panels), the 
back-scattering mechanism becomes 
more pronounced, creating the two 
red peaks. At small velocities, the 
results become similar to that for a 
static sphere. 



density and the expansion velocity. Note that, for small col- 
umn densities (upper panels), a certain fraction of photons 
were able to escape directly without scattering. These pho- 
tons form a very sharp peak at x = 0, and at small impact 
parameter, which has been removed from these plots for 
clarity. 

As the column density increases (top panels), the spec- 
trum becomes redder, and the small blue peak disappears, 
as the photons must scatter further from line centre in order 
to escape. The most significant change in the surface bright- 
ness profile occurs for small column densities, where some 
of the photons are able to escape directly, leaving fewer to 
scatter at large radii. 

As the expansion velocity increases (bottom panels), 
the back-scattering mechanism becomes more pronounced. 
Photons which scatter off the far side of the shell back 
through its interior are far enough from line centre to escape 
through the front side of the shell. This mechanism creates 
the two red peaks (the reddest peak comes from photons 
that backscatter more than once). At small velocities, we 
approach a profile similar to the static sphere. 

The surface brightness profiles are very flat when com- 
pared to the infall/outflow models of the gas in an NFW halo 
which we had considered before. Qua litati vely, they appear 
much fla tter t han the profiles of the iROSl emitters. Almost 
all of the lROSl emitters show a central peak. One way to pro- 
duce a central peak for a shell geometry is for the column 
density to be low enough for photons to be able to escape 
directly; the resolution of the instrument then broadens the 
delta function at y — into a central peak. This would limit 
the column density of the shell to A?m <2x 10 19 cm -2 . We 



conclude that the R08 emitters appear unlikely to be halo- 
scale expanding shells of Hi around a central Lya source. 



3.4 Summary of general trends of the Lya 
emission 

The radial (column) density distribution and the velocity 
field are the physical properties which most strongly affect 
the spectral distribution of the Lya emission in our mod- 
elling. For a spherical density distribution with neither out- 
flow nor infall, the distance between the peaks increases with 
increasing optical depth in the same way as the uniform 
static slab solution. Infall and outflow lead to a suppression 
of the red/blue peak and an increased shift of the opposite 
peak. The suppression increases with increasing velocity am- 
plitude and optical depth. 

The surface brightness profile depends strongly on the 
radial velocity profile. With increasing amplitude of the bulk 
motion, the diffusion in frequency space is accelerated and 
the emission becomes more centrally peaked. In our model 
with increasing amplitude towards decreasing radius (a = 
—0.5), the average bulk motions are larger and the effect is 
more pronounced. In the model with the expanding shell, 
the photons can travel unimpeded until they encounter the 
shell which leads to a rather flat surface brightness profile. 

The photons diffuse radially until the column density of 
neutral hydrogen density drops to values around 10 16 cm -2 
somewhat dependent on the spatial profile of neutral hydro- 
gen and the velocity field. The faint extended Lya emission 
should thus have a rather sharp edge which is defined by 
the surface inside which the gas is able to self-shield and the 
optical depth rises rapidly. 
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4 MODELLING THE RAUCH ET AL. 
EMITTERS 

4.1 Surface brightness profiles and spectral shapes 

We begi n with a qualitative summary of the properties of 
the IROSi emitters. Note, however, that due to the faintness 
of the sources the spectral and spatial profiles are rather 
noisy, making it difficult to identify signatures of inflow 
(moreprominent blue p eak) or outflow (more prominent red 
peakjj. As discussed in lROSl for 12 of the 27 spectra only a 
single emission peak is visible while six/three of the spectra 
show a weak secondary blue/red counter-peak. The remain- 
ing spectra are extended in frequency space without a clear 
peak structure. The widths of the spectral peaks ranges from 
~ 250 - 1000 km s"\ which corresponds to Aa; ~ 20 - 80 
for gas with a temperature of 10 4 K. The surface brightness 
profiles are predominantly centrally peaked with wings that 
often extend well beyond the Gaussian core of the PSF. This 
is particularly true of the brightest sources, while the fainter 
sources are more difficult to characterise due to the noise. 



4.2 Simultaneous modelling of DLA properties 

Modelling the radial distribution of neutral hydrogen in DM 
halos allows us to calculate the column density distribution 
along sightlines that intersect the halos. Comparison with 
the observed column distribution of DLAs provides a useful 
co nstrain t in this regard. We will use here the same notation 
as lBH09l . 

The column density distribution is defined such that the 
number of systems {d 2 M) intersected by a random line of 
sight between absorption distanc^fl X and X + dX, with Hi 
column density between JVhi and JVhi + diVm is, 



d 2 A/" = f(N,X) dX diVHi. 



(7) 



For our model we calculate this quantity using the Press- 
Schechter formalism. We need two ingredients. The first is 
um(M,X), the mass function of d ark matter halos, as cal- 
culated by ISheth fc Tormenl (|2002l ). The second ingredient 
is the column density of neutral hydrogen in a given halo 
(of mass M at absorption distance X), as a function of the 
(physical) impact parameter y, Nm(y\M,X). This is calcu- 
lated from the neutral hydrogen density as a function of 
radius. Given that Nm is a monotonically decreasing (and 
thus invertible) function of y, the region between Nm and 
TVhi + diVm is an annulus with cross-sectional area d(ny 2 ). 
Hence, we can write the column density distribution as, 



f(N,X) 



c 

Ho 



n M {M,X) 



d(7n/ 2 ) 



d7V H 



(Nm\M,X) 



dM. 



(8) 

The results are shown in Figure [7] The most important 
parameters (given that we integrate over all M) are the nor- 
malisation of the concentration parameter for the baryons 
Co, and the baryon fraction / e relative to the cosmic value. 

5 Recall also that this can be more complicated for the case of 
an outflowing shell, as there are then two red peaks. 

6 The absorption distance is defined by 

H 



dX 



H{z) 



(i + z y dz 



(6) 
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Figure 7. The column density distributio n f(N, X) for our 
model ; the black crosses show the data of IPr ochaska fc Wolfel 
(2009). The upper panel shows the effect of changing the con- 
centration parameter of the radial gas density profiles. The lower 
panel shows the effect of altering the baryon fraction f c (defined 
relative to the cosmic value). 



A good fit to the observed column density distribution is ob- 
tained for Co ~ 25.3 and f c ~ 0.2. The figure demonstrates 
the effect of altering these two parameters. Decreasing f c 
decreases the overall normalisation, while Co mostly affects 
the high column density end of the distribution. Changing 
v Ct o has a comparatively small effect. 

Note that our model somewhat overpredicts the num- 
ber of absorption systems with column densities below 
Am = 20.3 cm -2 , in the regime of Lyman Limit sys- 
tems (LLs). In this column density range the o bserved 
f{N,X) flattens significantly (|Q'Meara et al.ll2007l ). an ef- 
fect which can be attributed to the inability of (super) LLs 
to self-s hield completely against the met a-galactic UV back- 
ground (|Zheng fc Miralda-Escudell2002l ). Since we have not 
attempted to model the ionization of the gas and in par- 
ticular the self-shielding of the gas in any detail, it is not 
surprising that our model does not reproduce this. The UV 
background should, however, not effect our results for the 
population of predominantly neutral DLAs. Figure 5 shows 
that in the LLs regime the inferred size of the emission re- 
gion depends very weakly on column density. 
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Figure 8. The velocity width distribution /(ti w ,X) of the as- 
sociated low-ionization metal absorption of DLAs. The black 
crosses show the observational dat a compiled in Figure 10 of 
Wolf e. Gawiser. fc Prochaskal d2005h . The legend shows the pa- 
rameter v c fi, below which the baryonic fraction is assumed to 
be suppressed due to the effect of photo- heating and/or galactic 
winds. 



Similarly to IBH091 we can also calculate the predicted 
probability distribution of the velocity width u w of the as- 
sociated low ionization metal absorption of DLAs for our 
model as follows, 



10 



w,X) — — 



pOw I % (M) )hm (M, JQctdla (M, X)dM, 

(9) 

where p(v w \v c (M)) is the conditional proba bility distri- 
bution as discussed in Section 2 of IBH09I . The DLA 
cross-section is given by <jdla(M) = '^yhi.A, where 
Nth(vdla\M,X) = 10 20 - 3 cm' 2 . The result is shown 
in Figure El along with the observational data of 



IWolfe. Gawiser fc Prochaskal (|2005l ) 

As in BH09 our model fits the data well with values of 
v c .o in the range 50 — 70 km s _1 . We should point out that 
the probability distribution p(u w |fc) was originally derived 
from simulations which do not include the effect of galactic 
winds and where the distribution of gas in a given halo is 
somewhat different from what we have assumed here (see 
iBarnes fc HaehneltJ (|2009l ) for a more detailed discussion). 

We should also emphasize that, rather than using a 
simple powe r-law s caling for the absorption cross section of 
DLAs as in IBH091 we have used here a radial distribution 
of neutral hydrogen which is simultaneously consistent with 
the column density distribution of DLAs and the size dis- 
tribution of the Rauch et al. emitters, as we will see in the 
next section. 



4.3 The size distribution and luminosity function 

We will not attempt to fit the individual rath er noisy spec- 
tral and surface brightness profiles of the |RQ3 emitters here, 
but we will instead focus on the statistical properties of the 
population of emitters. 

We calculate the cross-section weighted size distribu- 
tion dN/dz (> r) expected from our model as follows. The 



Rauch et al 
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Figure 9. The cumulative size distribution dN/dz (> r) of the 
Lya emitters, compared with the observations of Rauch et. al 
(black curve). The coloured solid curves are for a = —0.5, while 
the dashed curves are for a = 1. The red and green lines are for 
the L oc Mm model, with the values of Lq as given in the legend. 
The blue curve assumes that r = r v , that is, assuming we can 
sec emission all the way to the virial radius. The lines have been 
normalised to dN/dz = 0.23 assuming a duty cycle of the Lya 
emission with /j, such that n^™ lttcrs = n^ los . The values of 



/d for each model (in the order they appear in the legend) are 
/ d = 1,0.2,0.055 for a = -0.5 and f d = 0.28,0.07,0.055 for 
a = 1. The models with the two values in bold correspond best 
to the data. 



observations of iRauch et al.l (|2008l ) achieved a 1 a surface 
brightness detection limit of So = 10~ 19 erg s -1 . We calcu- 
lated the expected observed size of our model emitters by 
determining the radius r (or equivalently, impact parame- 
ter) at which the surface brightness drops below the iROSl 
limit, S(r) — So- 

This procedure gives the radius of the emitter as a func- 
tion of the mass of the halo. To do this, we need to specify 
the intrinsic Lya luminosity L^ ya as a function of mass. 
Similarly to[BH09, we assume that the luminosity is pro- 
portional to the total mass of neutral hydrogen, 



= L 



erg s 



(10) 



,2.4 x 1O 9 M 

This means that the luminosity is subject to the same expo- 
nential suppression as the (neutral) gas content of DM halos 
for small circular velocities. 

We now calculate the size distribution in the form of the 
inferred cumulative incid ence rate dN/dz (> r) and compare 
it with the data of lR08l in Figure [9] The curves have been 
normalised to dTV /dz — 0.23 by assuming the emission oc- 
curs with a duty cycle / d , n e ^ ttela = ,f d n\f os . The values 
of fd for each model are given in the caption to the figure. 
The solid curves are for a = —0.5, while the dashed curves 
are for a = 1. The red and green curves assume L oc Mm, 

7 The constants of proportionality are chosen to that without the 
suppression below ti c ,p and with f e = 0.2, Lq defined below has 
the same value as in iBHOSt 
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Figure 10. The cumulative luminosity function n(> LLya) of 
our model, along with the Rauch et al. data (black solid curve 
with errors indicated by the grey dot-dashed curves). The same 
models are shown as in Figure [9] except that the duty cycle is 
chosen to match the luminosity function data. The parameters 
for each model (in the order they appear in the legend) are (with 
negligible dependence on a): fd = 0.28,0.07. 



with the values of Lo as given in the legend. The blue curve 
shows the size distribution assuming that r — r v , that is, 
assuming that the emission is detected all the way to the 
virial radius. 

As we have defined the size of the emitters at a fixed 
flux level So, the observed size of the halos becomes larger as 
Lo increases. For a = —0.5 the region with emission above 
this flux level is smaller than in the model with a — 1 halos, 
due to the more centrally concentrated surface brightness 
profile. The size distribution flattens at small r due to the 
exponential suppression of the luminosity at low v c . 

Figure [9] shows that we can find values of Lo that fit the 
observed cumulative size distribution well, where the value 
of Lo depends on a. The required duty cycle is ~ 20 — 28%, 
and is fairly insensitive to a. The value of Lo for a = —0.5 
is rather high, due to the very peaked surface brightness 
profiles in this case. 

Figure [10] shows the corre spond ing luminosity distribu- 
tion n(> Z/Lya), along with the IROSi data (black solid curve). 
Note that the luminosity Z/L y a predicted by our model is 
calculated by integrating the surface brightness inside the 
radius r where the surface brightness is above the observa- 
tional limit. Photons that are scattered to radii where the 
surface brightness falls below this limit are lost in the noise; 
the observed luminosity is thus always less than the intrin- 
sic luminosity. Note further that the observed luminosities 
have not been corrected for "slit losses". The actual Lya 
luminosities may thus be a factor two or more larger. Note 
finally, that the intrinsic Lycv emission could be significantly 
larger due to absorption by dust. 

The figure shows the same models as in Figure [5] except 
that the duty cycle is chosen to fit the luminosity function 
data. The solid and dashed curves are very close, meaning 



that there is practically no dependence on a. The corre- 
sponding values of fa are given in the caption to the figure. 

The observed lumininosity function is reproduced rea- 
sonable well by our model if we assume Lo = 10 42 erg s — 1 , 
the value required to match the observed size distribution 
for our model with a = 1 This model predicts a somewhat 
steeper fain t end slope than appears to be observed, but as 
discussed in R08 it is uncertain whether the turn-over in the 
observed luminosity function is real or due to incomplete 
identification of emitters close to the detection threshold. 
Should the turnover consolidate with deeper data then this 
may suggest a somewhat sharper cut-off of the efficiency for 
Lya emis sion in shallow potential wells as we had e.g. as- 
sumed in BH09. The required duty cycle is fd = 0.28 and 
thus comparable to the value required to fit the observed 
size distribution. With a = I, L = f0 42 erg s _1 , f d ~ 0.28 
our model therefore f its bo th the observed size and luminos- 
ity distribution of the lR08l emitters. As previously discussed, 
the model with a — —0.5 has strongly centrally peaked emis- 
sion. As is apparent from Figure [T0]this leads to a mismatch 
with the observed luminosity function if we fix Lo = fO 43 ' 3 
to fit the observed size distribution for this value of a. Even 
for a significantly smaller duty cycle of fd = 0.07 we cannot 
match the shape of the observed luminosity function, which 
in this case is significantly steeper than predicted by the 
model over the full range of luminosities. 

We therefore conclude that we can successfully repro- 
duce both the absorp tion properties of DLAs and the Lya 
emission data of lRO&l using a self-consistent model with cen- 
trally concentrated, star-formation powered Lya production 
with a duty cycle of ~ 25%, coupled with radiative trans- 
fer effects that set the observed size of the emitters if the 
velocity field of the gas facilitates the scattering of the pho- 
tons to large radii with moderate central bulk velocities as 
in the model with a = f . Our model's success in simultane- 
ously reproducing the absorption p roperties of DLAs s ignif- 
icantly strengthens the assertion of iRauch et al.l (2008) that 
the faint emitters are in fact the host galaxies of DLAs. 

Finally we calculate d 2 J\f /dX/d log 10 M, the contribu- 
tion to the incidence rate from halos of different masses. The 
result is shown in Figure [TTl The majority of DLAs/emitters 
have masses in the range 10 9 ' 5 — 10 12 Mq for all models. This 
is a similar range of masses to that fou n d in the numeri- 
cal simulations o f Nagamine et al. (2007); Razo umov et al.l 
(|2008t ); iPontzen etail (120081 '); iT^scari et all (|2009l ). We re- 
fer the reader to |BH09| for a more detailed comparison of 
our prediction of d 2 A/"/dX/dlo g ln M with th a t of n umeri- 
cal simulations. The success of IPontzen et al] |20 06f) in re- 
producing the observed metalicity distribution is additional 
evidence in favour of this mass range, as are the observe d 
spatial correlation of DLAs and LBGs (|Cooke et al.l l2006). 

4.4 Differences to our previous modelling 

In lBHQgl . we presented modelling in which we assumed sim- 
ple scaling relations for the sizes and luminosities of the 
Rauch et al. emitters with the masses/circular velocities 
of their DM host haloes. The scaling relation for the sizes 
was motivated by the modelling of the velocity width dis- 
tribution of the associated low ionization m e tal ab sorption 
of DLAS by lHaehnelt. Steinmetz. fc Rauch! (|f 99&f ). No at- 
tempt to model the spatial and spectral distribution of the 
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Figure 11. The contribution of different mass ranges to the 
incidence rate of DLAs/emitters (d 2 A/"/dX/d log 10 M). The same 
models are shown as in Figure [9] The majority of DLAs/emitters 
have m asses in the range 10 10 — f0 12 MQ. We refer the reader to 
IBH09I for a comparison of our prediction of d 2 Af /dX/d log 10 M 
to that of numerical simulations. 

Lya emission was made in |BH09| . It is thus gratifying that 
our modeling here with Lya radiative transfer through in- 
flowing/outflowing distributions of neutral hydrogen in the 
emitters for which the radial profile is constrained by the 
observed column density ditribution is able to still fit the 
luminosity function and size distribution of the emitters as 
well as the velocity width distribution of DLAs for a specific 
choice of the scaling of the inflow/outflow velocity with ra- 
dius. This is a non-trivial result and means that the radiative 
transfer calculations as constrained by the observed proper- 
ties of the Rauch et al. emitters and DLAs produce a scaling 
relation of the sizes of emitters with the masses/circular ve- 
loc ities of their DM host haloes very simular to that assumed 
bv lBH09l . 



4.5 Limitations of the current modelling 

Our modelling has a number of limitations and (over-) sim- 
plifications that we will discuss in this section. Firstly, we 
have assumed spherical symmetry throughout. Deviations 
from a uniform spherical configuration will most likely make 
escape easier for Lya photons in certain directions. Simi- 
larly, our simple velocity field is certainly an oversimplifica- 
tion. More complex velocity fields may bring photons pass- 
ing through gas in the outer regions of the halo back toward 
line centre. This may also increase the spatial extent of the 
emitter. 

We have also ignored radiative transfer outside the 
virial radius of the halo. This effectively assumes not only 
that scattering in the outer parts of DM halos is negligible 
but that the same is true for scattering by the neutral hy- 
drogen contained in the IGM. Whether these are reasonable 
assumptions will depend both on the poorly constrained dis- 
tribution of neutral hydrogen in the outer parts of low-mass 



DM halos and on the bul k motion of th e scattering gas rel- 
ative to the IGM (see e.g. ISantosl lj2004) h 

We have furthermore not modelled the effect of ion- 
ising radiation, either from external or internal sources of 
ionising photons, except by reducing the amount of neutral 
hydrogen within the virial radius of the DM halos below its 
cosmic value, ostensibly due to the effects of the UV back- 
ground. More sophisticated modelling should take into ac- 
count the self-shielding of the gas against the meta-galactic 
UV background self-consistently. We have checked here that 
the radius at which the gas is able to maintain neutrality by 
self-shielding is > r om i ss i on , the size of the emission region 
at the flux level of the Rauch et al. emitters as defined in 
Section 14.31 The meta-galactic UV background should thus 
have little effect on the size of the Lya emission region pre- 
dicted by our modelling. Note, however, that this will not 
be the case anymore at fainter flux levels as the neutral hy- 
drogen density /column density will drop sharply outside the 
region able to self-shield. 

Internal sources of ionising stars such as stars would 
reduce the neutral fraction of the gas, and in an inhomo- 
geneous manner. Neither stars nor the UV background are 
expected to significantly ionise the bulk of the neutral gas 
in the proto-galaxies studied here — these are the sites of 
DLAs after all. Spatially extended star formation is likely 
to increase the production of Lya at large radii, possibly 
increasing the observed size of the emitters. 

Finally, we have also ignored the effects of dust. Our 
assumption of a negligible effect of dust may not be unrea- 
sonable as DLAs are known to have a rath er small dust 
content IjWolfe. Gawiser. fc Prochaskal l2005t) and the es- 
cape fraction for Lya emission appears to increase strongly 
towar ds unity with d ecreasing stellar mass of the emit- 
ters (lOno et al.l 120091). The recent modelling of Lya ra- 
diativ e transfer by Laursen. Sommer-Larsen. fc AndersenI 
(2009) including dust also suggest that for the majority of 
DLAs, with their rather low metallicity which are hosted by 
rather low mass halos, the effect of dust is probably not too 
important. Note however that there are coun ter-examples. 
As shown bv lAtek. Schaerer. fc Ku nth (2009), models with 
almost complete absorption of dust can explain the fact that 
the expected Lya emission from the star formation in IZw 
18 appears to be entirely absorbed despite a rather mod- 
est reddening. As discussed in section 3. IT, the presence of 
dust may actually be (partially) responsible for the smaller 
than unity duty cycle which we require to fit the luminosity 
function. Note also that the effect of dust shou ld depend sen- 
sitive ly on the dumpiness of the neutral gas (|Hansen fc OhI 
l2006h . 



5 CONCLUSIONS 

We have used a Monte-Carlo radiative transfer code to 
model the spatial and frequency distribution of the Lya 
emission due to star formation in (proto-)galaxies at the 
centre of DM halos with masses of I0 9 ' 5 to I0 12 M Q for a 
range of assumptions for the spatial distribution and the 
dynamical state of neutral hydrogen. DM halos in this mass 
range had been previously identified as the likely hosts of 
DLAs and the recently detected population of faint spatially 
extended Lya emitters. Our main results are the following: 
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• As previously found by other authors, the spectral 
shape of the Lya emission from star formation in galaxies 
for which the dynamics of neutral hydrogen is dominated 
by infall/outflow is characterized by a strong blue/red peak, 
occasionally accompanied by a weaker red/blue peak. The 
spectral shape is very sensitive to the spatial distribution, 
velocity structure and (to a smaller extent) the tempera- 
ture structure of the gas. The larger central column densi- 
ties in the more massive galaxies/ halos make escape of the 
Lya photons more difficult. The photons must then scat- 
ter further in frequency and space and thus emerge with 
a larger frequency shift and the emission extends to larger 
radii. Larger bulk motions lead to more energy per scatter- 
ing being transferred between the photons and the gas. This 
results in a larger frequency shift of the dominant spectral 
peak and a larger contrast between the strong and weak 
spectral peak. 

• The surface brightness profiles for photons emitted at 
the centre of the halos show a central peak with wings ex- 
tending as far as our assumed neutral hydrogen distribution 
as long as the column density of neutral hydrogen exceeds 
about ~ 10 16 cm~ 2 . The spatial profile of the emission is 
likewise sensitive to the spatial distribution, velocity and 
temperature structure of the gas. The spatial distribution is 
significantly more centrally peaked when the amplitude of 
the bulk motions increase toward the centre of the halo. 

• Expanding shells of neutral hydrogen similar to those 
invoked to explain the Lya emission from LBGs produce 
spectra with one or more prominent red peaks. The surface 
brightness profiles are very flat, remaining essentially con- 
stant for 75% of the radius of the shell . This appears at odds 
with the observed profiles of the lROSI emitters, almost all of 
which show a central peak. 

• Our modelling simultaneously reproduces the column 
density distribution of the neutral hydrogen and the veloc- 
ity width of the associated low ionisation metal absorption 
of DLAs, as well as the size distribution and the luminosity 
function of the Rauch et al. emitters if we assume i) that 
absorbers and emitters are hosted by DM halos that retain 
about 20% of the cosmic baryon fraction in the form of neu- 
tral hydrogen, with a spatial distribution which follows a 
NFW profile with concentration parameter ~ 7 times larger 
than that of the dark matter, ii) that absorbers and emitters 
are hosted by DM halos with virial velocities > 50 km s" 1 , 
and iii) that the central Lya emission due to star formation 
has a duty cycle of ~ 25% and the luminosity is proportional 
to the mass of neutral hydrogen in the DM halos. 

• The DM halos that contribute most to the incidence 
rate of DLAs have masses in the range M ~ 10 9,5 - 1O 12 M 
and virial velocities in the range of 35 to 230 km s" 1 . The 
lower cut-off is mainly determined by the rather sharp de- 
crease in the velocity width distribution of the associated 
low-ionization metal absorption in DLAs at velocity width 
< 30 km s _1 but may also be reflected in the turn-over of 
the Lya luminosity function at the faintest fluxes. The DM 
host halo masses are significantly smaller than those inferred 
for L» LBGs, consistent with the much higher space density 
of the faint emitters. 

The success of our detailed Lya radiative transfer mod- 
elling in explaining the observed properties of both DLAs 
and the faint Rauch et al emitters with a consistent set 



of assumptions further strengthens the suggestion that the 
faint emitters are indeed the long-searched for host galax- 
ies of DLA/LLs. Together with our modelling, the observed 
properties of the faint emitters should thus provide robust 
estimates not only of the space density, Lya luminosity and 
extent of the gas distribution but also of the masses of the 
DM halos, the duty cycle of star formation and the spa- 
tial profile and kinematics of the gas distribution of a sta- 
tistically representative sample of DLA host galaxies. The 
current ultra-deep spectroscopic surveys in the HUDF and 
HDF should soon provide important additional information 
on the stellar content and possibly also dust content of these 
objects which will allow to further test the nature of what 
are almost certainly the building blocks of typical present- 
day galaxies like our own. 
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APPENDIX A: Lya RADIATIVE TRANSFER 
ALGORITHM 

Al Structure of the Monte-Carlo Code 

The details of each step in our Monte Carlo code for Lya 
RT are outlined below — throughout, Ri (for i = 1,2,3,...) 



denotes a random number generated uniformly between 
and 1. 

1. ) We begin by specifying, as a function of position, r 

- nm(r), the number density of Hi, 

- e(r), the Lya emissivity (in photons/s/cm 3 ), 

- Vb(r), the bulk velocity field of Hi, and 

- T(r), the temperature of Hi. 

For the modelling in this paper, these quantities are dis- 
cussed in Section \3. II 

2. ) We generate a photon at an initial position r; accord- 
ing to the emissivity. We then choose the photon's initial 
direction ft from an isotropic distribution, 



n = (sin 9 cos <f>, sin 9 sin <f>, cos 9) 



(Al) 



where 



= cos~ 1 (2i?i - 1) 
' = 2nR 2 



(the polar angle), (A2) 
(the azimuthal angle). (A3) 



We generate the photon's initial frequency as follows. The 
emission (and absorption) profile in the rest frame of the 
emitting atom is assumed to be Lorentzian, 



0i» = 



A^l/2tt 



(y - v y + (A^l/2) 2 



(A4) 



where v$ — 2.47 x 10 15 Hz is the central Lya frequency and 
Az/l = 9.936 x 10 Hz is the natural line width. We take into 
account the thermal Doppler broadening with an Maxwell- 
Boltzmann (MB) distribution of the velocity of the scatter- 
ing atoms. The line profile can then be written as an average 
of the Lorentz profile over the atoms' velocities, 



oo 

f (l 



1 



exp 



H(a,x), 



(A5) 
(A6) 



where the argument of </>l takes into account the Doppler 
shift in the frequency of the photon as seen by an atom 
with velocity component v z in the direction of the photon's 
motion. The other quantities are 



Vth = 



(m 1/2 = 12.85 km s- f^Y 



\ m j 
Av h 



V 10 4 K 



2Ai/d 

'Vth 



4.693 x 10" 



/ T 



V10 4 K 



Ai/d = ^— ) VQ, 
v - v Q 



(A7) 

(A8) 
(A9) 
(A10) 



the thermal velocity dispersion, relative line width, Doppler 
frequency width and frequency displacement in units of 
Also, respectively. We will use x as our frequency variable. 
We have also made use of the Voigt function, 



H(a, x) 



(x — y) 2 + a- 



Ay. 



(All) 
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The Voigt function can be approximated as a Doppler core 
and Lorentz wings, 



H(a, x) 



e x , \x\ < x c (core), 
7=^7, \x\>x c (wing). 



(A12) 



The transition between the two occurs at a 'critical' fre- 
quency x c , defined as the solution to e~ x " = a/^nXc- For 
T = 10 4 K, we have x c ~ 3.255. W e will use the app roxima- 
tion to the Voigt function given in iTasitsiomil (|2006l ). 

Equation (IA5|) assumes that the "laboratory frame" (in 
which we are measuring the photon's frequency) is the same 
as the fluid rest frame. In order to take account of the bulk 
fluid velocity, we replace the frequency in the lab frame, x, 
with the frequency as seen by the fluid, 



Vb ■ n 

Vth 



(A13) 



In practice, the photons are very likely to be emitted close 
to line centre, and because the photons suffers a large num- 
ber of resonant scatterings, any "m emory" of the initial fre- 
quency is quickly erased (jDHSOd ). Thus, we will usually 
inject all photons at line centre in the fluid frame (x — 0). 

3.) The distance travelled depends on the optical depth. 
The probability that the photon propagates a physical dis- 
tance corresponding to an optical depth between r and r+dr 
is e~ T dr. We choose an optical depth for the photon from 
this distribution, 



HR 3 



(A14) 



To find the physical distance travelled, we perform the 
following line integral, solving for Sf, 



c(r(s)) nm(r(s)) ds, 



(A15) 



where r(s) = Vi + hs is the path travelled by the photon, 
and a x is the scattering cross-section of Lya photons, 



fi2itcr e (f)(x) 



(A16) 



and /12 = 0.4167 is the Lya oscillator line strength while 
r e = 2.81794 x 10 -13 cm is the classical electron radius. 
It is worth noting how the integrand depends on position. 
The most obvious dependence is that of a spatially varying 
density nm(r). The dependence of o v is two-fold. Firstly, 
in case of a spatially-varying bulk velocity, a v [cj>{x{vb{r)))). 
Secondly, if T depends on position, then so will both x and 
a in (f>{x) via Ai^d and vth- 

Once Sf is found, the position of the next scattering is 
r = Ti + nsf. If this is outside the Hi region, then the 
algorithm is terminated. 

4.) Next we choose the velocity of the scattering atom. 
Naively, one might think that this step involves generat- 
ing the three components of the atom's velocity from the 
Maxwell-Boltzmann distribution. However, we are choosing 
the velocity of an atom given that it scatters a photon with 
frequency x. We can therefore divide the velocity of the atom 
into one component parallel («,,) and two components per- 
pendicular to the direction of motion of the photon («xi, 
v±2)- The two perpendicular components do not alter the 
frequency of the photon as seen by the atom, and are thus 
chosen from a Maxwell-Boltzmann distribution. From Equa- 
tions (|A5[) - (|A11|) ^ can be seen that the probability that a 
photon with frequency x scatters off an atom with velocity 



(along the direction of propagation of the photon) between 
v„ and v„ + dv„ is, 



/(«„) du„ = - 



tt (x — M„) 2 + a 2 H(a,x] 



dun 



(A17) 



where u„ = v„/v t ^. A scheme for generating ran- 
dom numbers from t his distribution is given in 
I Zheng fc Miralda-Escudl <|2002h ■ 

5.) Now that we have calculated v a , the velocity of the 
atom that scatters the photon, we perform a Lorentz trans- 
form into the rest frame of the atom. In this frame, we as- 
sume that the frequency of the scattered photon differs from 
the f requency of the incident p hoton only by the recoil ef- 
fect (|Rvbicki fc Lightman|[l979l . pg. 196). We choose a new 
direction for the photon. Our code can incorporate either 
an isotropic or a dipole distribution for the direction of the 
re-emitted photon. The results are insensitive to the choice. 
We then transform back into the laboratory frame. At speeds 
much less than c we have, for an initial photon direction n, 
and having chosen a final photon direction n/, that the final 
frequency of the photon is, 

n-V a fif-Va hAu D „ . . . . 

x f = x + ^ + — — (n • n f - 1). A18 

«th «th 2fcJ 

In ge neral, the final term, known as the recoil term (jAdamsl 
Il97lf ). is negligible for the modelling in this paper. 

We now return to step (iii) and repeat it until the photon 
escapes the Hi region. 

Once the photon escapes the region, its properties (fre- 
quency, angle of escape etc.) are recorded. We then return 



to step (ii) and generate another photon. 



Our code incorporates the presence of a cosmic abun- 
dance of deuterium, following the method presented in 
DHS06. Adding deuterium has only a minimal effect. 



A2 Accelerating The Code 

Monte Carlo Lya RT codes like the one used here can be sig- 
nificantly accelerated by skipping scatterings in the core of 
the line profile. We define the core-to-wing transition to oc- 
cur at a critical frequency Xctit, which is not the same as x c . 
Whenever a photon is in the core \x\ < x CI u, we force it into 
the wing by choosing the scattering ato m's velocity to be 
large. We follow the method presented in lAhn. Lee, fc Leg 
(2002), DHS06, choosing the perpendicular components of 
the atom's velocity such that Uj_i + u \i 1? x 1mt- 

We choose the value of x cr it by requiring that, for a 
uniform sphere, less than a fraction / of the photons that 
emerge have \x\ < Kcrit- Using the analy tic solution for the 
emergent spectrum of a sphere (|DHS06t Equation (9)), this 
gives, 



/54 ♦ k-i / 
— aro tanh / 



(A19) 



where to is the line-centre optical depth from the centre to 
the edge. We find that setting / = 0.01 up to max(i cr i t ) = 
3 (~ x c ) gives an acceptable compromise between speed and 
accuracy. This procedure can accelerate the code by several 
orders of magnitude, e.g. for a static sphere with To = 10 7 , 
the code is accelerated by a factor of more than a thousand. 
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A3 Spherically Symmetric Shells 

In cases where nm, T and/or Vh have a complicated depen- 
dence on r, solving for s/ in Equation (|A15[l can be compu- 
tationally ex pensive. In spherical symmetry, we use a similar 
approach to iDHSOq , dividing the sphere into shells of uni- 
form density. Within each shell, the integral (|A15[) becomes 
trivial, r = s/<7 x nHi(r s heil)- If the edge of a shell is reached, 
the optical depth to the edge of the shell is subtracted from 
r and a new Sf is calculated using the new value of r s h n. 
We space the shells so that each shell has equal column den- 
sity. We choose the number of shells to keep the frequency 
dispersion within each shell small compared to the change 
in frequency as the photon crosses a shell; 1000 shells are 
usually sufficient. 



APPENDIX B: THE EFFECT OF 
TEMPERATURE AND VELOCITY 
AMPLITUDE 

In this section, we will consider the effect of altering the gas 
temperature and amplitude of the bulk velocity. For this 
section, the fiducial model is different to that in Section 
13X41 (z,M„,co,/ e ,«a m p,T) = (3, 1O 11 M ,3.5,1,1) C ,1O 4 K). 

Bl The effect of the velocity profile 

In this section, we consider the effect of decreasing v a mp- 
The results are shown in Figure iBll Setting u amp = would 
result in J(x) being symmetric about x = (ignoring the 
small effects of deuterium and recoil). This tendency is seen 
clearly in both of the left plots. These plots also show that 
photons emerge bluer as the velocity of infall increases. Pho- 
tons gain energy from head-on collisions with atoms, and the 
more energetic the atoms, the more energy is transferred be- 
tween photons and gas. 

The surface brightness plots show that the emission be- 
comes more extended as the velocity is reduced. A rather flat 
surface brightness profile is characteristic of a uniform, static 
Hi sphere. This is because the bulk velocity can give pho- 
tons a "free ride" through the halo, Doppler shifting their 
frequency away from line centre in the fluid frame without 
changing x in the laboratory frame. In the a = —0.5 case, 
the peak surface brightness is reduced by almost 2 orders 
of magnitude by this effect. It is also worth noting that the 
smaller « a m P is, the less the dependence of the profile on a. 

B2 Temperature 

The effect of altering the temperature is shown in Figure 
IB2I Note that the spectrum is given as a function of ve- 
locity (v = —x Vth), as x itself is temperature dependent 
(oc 1/VT). For the a = -0.5 case, there is no dependence 
on the temperature because the photon begins its flight with 
a frequency (in the fluid frame) in the scattering wings, and 
is unlikely to return to the Doppler core. Remember that 
the scattering cross-section does not depend on T in the 
Lorentz wing of the profile. In the a = 1 case, the effect of 
temperature is minimal. The red peak disappears at lower 
temperatures due to the corresponding increase in the scat- 
tering cross-section, a x oc 1/Auo oc 1/vT, when x is small. 



This paper has been typeset from a 'TpX/ I^TpX file prepared 
by the author. 
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Figure Bl. Spectra (left pan- 
els) and surface brightness pro- 
files (right panels) for DLAs in ha- 
los with Damp as given in the leg- 
end in units of the virial veloc- 
ity. The top panels are for a = 
— 0.5, for which L^ ya = 1.1 X 
10 44 erg s — 1 . The bottom panels 
are for a = 1, for which Li^ ya = 
5.5 X 10 42 erg s _1 . As the velocity 
decreases, the spectrum begins to 
resemble the double-peaked static 
sphere profile. The surface bright- 
ness profile flattens as the veloc- 
ity decreases — where the bulk 
velocity is low, the photons must 
random walk out of the Doppler 
core by scattering; they will not 
be given a "free ride" by the fluid 
flow. 
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Figure B2. Spectra and surface 
brightness profiles as a function of 
velocity for DLAs in haloes with 
temperature To as given in the leg- 
end. The top panels are for a = 
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500 




-0.5, for which Llv 



1.1 X 



10 44 erg s — 1 . The bottom panels 
are for a = 1, for which Li^ ya = 
5.5 X 10 42 erg s _1 . For the a = 
—0.5 case, there is no dependence 
on the temperature because the 
photon begins its flight with fre- 
quency (in the fluid frame) in the 
scattering wings, and is unlikely to 
return to the Doppler core. In the 
a = 1 case, the effect of tempera- 
ture is minimal. The red peak dis- 
appears at lower temperatures due 
to the corresponding increase in 
the scattering cross-section, cr x oc 
l/Au-£, oc l/y/T, when x is small. 
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